1 About

This document contains the analyses of the manuscript entitled “Regional heterothermy is not associated with active heat dissipation by the horns in the Rhinoceros Beetle (Megasoma gyas)” by Danilo Giacometti, Luiz Henrique Lima Silva, Guilherme Gomes, José Eduardo de Carvalho, and Alexandre V. Palaoro.


2 Abstract

To maintain thermal balance, animals typically rely on the control of physiological and behavioral processes. Some animals, however, possess exaggerated morphological structures that aid in thermoregulation, particularly in the dissipation of excess heat when body temperatures rise. Although widespread in the animal kingdom, the role of animal weapons—exaggerated morphological structures with multiple characteristics—has received little attention in thermoregulation studies. Here, we investigated whether the horns of the Rhinoceros Beetle (Megasoma gyas) acted as a thermal window. We heated live and dead beetles to 30 ºC and allowed them to cool to 20 ºC, while measuring surface temperature changes in four body regions: the cephalic horn, the thoracic horn, the elytra, and the abdomen. If the horn actively dissipated heat, the cephalic horn would show the lowest cooling rate among body regions. Contrary to this expectation, we found that the cephalic horn had the highest cooling rate, followed by the abdomen, thoracic horn, and elytra, respectively. This suggests that the horns are not used actively for heat dissipation in M. gyas. Live specimens had lower cooling rates than control specimens across all body regions, which could be an indication of hemolymph flow modulating heat dissipation. However, whether beetles show circulatory responses to temperature is currently unknown. The low cooling rate of the elytra can be explained by the high density of flight muscles in the thorax. Apart from their role in heat generation, flight muscles may aid in heat dissipation by pumping hemolymph across tagmata or through the low-insulated cuticle to prevent thoracic overheating. Our study also demonstrates that invertebrates show regional heterothermy even in instances unrelated to exercise or stress. We propose that regional heterothermy results from both active (control of hemolymph flow) and passive (heat dissipation through poorly insulated structures) processes within individual.

3 Packages

library(scales)
library(performance)
library(lme4)
library(tidyverse)
library(forecast)
library(ggplot2)
library(kableExtra)
library(cowplot)
library(pander)
library(png)
library(patchwork)
library(sjPlot)

4 Functions

4.1 Standardised theme for plotting

ggtheme <- function(base_size=12, base_line=0.3) {
  theme(
    
    text =        element_text(size=base_size),
    line =        element_line(size=base_line, linetype="solid"),
    
    axis.text.x = element_text(size=base_size*0.8, colour='black',  hjust=0.5, vjust=1, angle=0),
    axis.text.y = element_text(size=base_size*0.8, colour='black', hjust=1, vjust=0.5, angle=0),
    axis.line.x = element_line(size=base_line),
    axis.line.y = element_line(size=base_line),
    
    axis.title.x =  element_text(size = base_size, vjust = 1, margin=unit(c(3,0,0,0),"mm")),
    axis.title.y =  element_text(size = base_size, angle = 90, vjust = 0.5, margin=unit(c(0,3,0,0),"mm")),
    axis.ticks = element_line(size=base_line),
    axis.ticks.length = unit(0.3, "lines"),
    
    panel.background = element_blank(),
    panel.grid.minor = element_blank(),
    panel.grid.major = element_blank(),
    panel.border     = element_rect(fill=NA, colour=NA, size=base_line),
    panel.spacing    = unit(1, "lines"),
    
    legend.background = element_rect(fill="transparent", colour="transparent"),
    legend.key =       element_rect(fill="transparent", colour="transparent"),
    legend.text=       element_text(size=base_size),
    
    strip.background =  element_blank(),
    strip.text.x =      element_text(size = base_size * 0.8),
    strip.text.y =      element_text(size = base_size * 0.8, angle = -90),
    strip.switch.pad.grid = unit(0, "mm"),
    strip.switch.pad.wrap = unit(0, "mm"),
    
    plot.margin =       unit(c(0.5, 0.5, 0.5, 0.5), "cm"),
    plot.title =        element_text(size = base_size * 1.2),
    plot.background =   element_rect(colour = "transparent", fill="transparent", size=base_line)
  )}

**Figure 1**. Body parts of *Megasoma gyas* considered in the current study. From each body part, we obtained a measure of surface temperature over five minutes after animals were passively heated to 30 ºC and allowed to cool down to room temperature (~ 20 ºC). We used the change in surface temperature over time to obtain a measure of cooling rate for each body part.

Figure 1. Body parts of Megasoma gyas considered in the current study. From each body part, we obtained a measure of surface temperature over five minutes after animals were passively heated to 30 ºC and allowed to cool down to room temperature (~ 20 ºC). We used the change in surface temperature over time to obtain a measure of cooling rate for each body part.


**Figure 2**. Temporal changes in heat dissipation across body parts in *Megasoma gyas*.We used a water bath to heat beetles to 30 ºC and allowed them to dissipate heat to their surroundings (kept at 20 ºC) for 5 min. Each letter depicts a 1-min increment from the (A) start (t = 0 min) until the (F) completion (t = 5 min) of the experiment. The cephalic horn cooled faster than the other body parts considered in this study.

Figure 2. Temporal changes in heat dissipation across body parts in Megasoma gyas.We used a water bath to heat beetles to 30 ºC and allowed them to dissipate heat to their surroundings (kept at 20 ºC) for 5 min. Each letter depicts a 1-min increment from the (A) start (t = 0 min) until the (F) completion (t = 5 min) of the experiment. The cephalic horn cooled faster than the other body parts considered in this study.


5 Data set processing

beetle<-read.csv(file = "cool-heat.csv", h = T, sep = ';')

beetle$invert<-beetle$cool.rate*-1 # obtaining positive cooling rate values
beetle2<-beetle[!(beetle$exp == "active.heating"),] # removing active heating experiments from data set

5.1 Creating subsets

beetle.control<-subset(beetle2,exp=="control")

beetle.passive<-subset(beetle2,exp=="passive.heating")

5.1.1 Control data set

kable(beetle.control) %>% kable_styling() %>% scroll_box(width = "100%", height = "500px")
ID exp body.part cool.rate mass mass.cat horn.ratio invert
49 C1 control cephal.horn -1.0610 22.86 high 0 1.0610
50 C1 control thorax.horn -0.3569 22.86 high 0 0.3569
51 C1 control elytra 0.0000 22.86 high 0 0.0000
52 C1 control abdomen -0.3344 22.86 high 0 0.3344
53 C2 control cephal.horn -0.6196 22.86 high 0 0.6196
54 C2 control thorax.horn -0.7033 22.86 high 0 0.7033
55 C2 control elytra -0.3986 22.86 high 0 0.3986
56 C2 control abdomen -0.2050 22.86 high 0 0.2050
57 C3 control cephal.horn -0.5532 22.86 high 0 0.5532
58 C3 control thorax.horn -0.9675 22.86 high 0 0.9675
59 C3 control elytra -0.9305 22.86 high 0 0.9305
60 C3 control abdomen -1.0210 22.86 high 0 1.0210

5.1.2 Passive heating data set

kable(beetle.passive) %>% kable_styling() %>% scroll_box(width = "100%", height = "500px")
ID exp body.part cool.rate mass mass.cat horn.ratio invert
1 passive.heating cephal.horn -0.3747 23.79 high 0.9540217 0.3747
1 passive.heating thorax.horn -0.1000 23.79 high 0.9540217 0.1000
1 passive.heating elytra -0.2329 23.79 high 0.9540217 0.2329
1 passive.heating abdomen -0.0839 23.79 high 0.9540217 0.0839
2 passive.heating cephal.horn -0.1852 22.90 high 0.6802248 0.1852
2 passive.heating thorax.horn -0.0521 22.90 high 0.6802248 0.0521
2 passive.heating elytra -0.0243 22.90 high 0.6802248 0.0243
2 passive.heating abdomen -0.1726 22.90 high 0.6802248 0.1726
3 passive.heating cephal.horn -0.3769 19.13 low 0.7269632 0.3769
3 passive.heating thorax.horn -0.0227 19.13 low 0.7269632 0.0227
3 passive.heating elytra -0.0412 19.13 low 0.7269632 0.0412
3 passive.heating abdomen -0.1996 19.13 low 0.7269632 0.1996
4 passive.heating cephal.horn -0.3623 23.32 high 0.6548562 0.3623
4 passive.heating thorax.horn -0.0675 23.32 high 0.6548562 0.0675
4 passive.heating elytra -0.0229 23.32 high 0.6548562 0.0229
4 passive.heating abdomen -0.3104 23.32 high 0.6548562 0.3104
5 passive.heating cephal.horn -1.2905 13.74 low 0.4223241 1.2905
5 passive.heating thorax.horn -0.1514 13.74 low 0.4223241 0.1514
5 passive.heating elytra -0.0918 13.74 low 0.4223241 0.0918
5 passive.heating abdomen -0.1208 13.74 low 0.4223241 0.1208
6 passive.heating cephal.horn -0.2828 22.86 high 0.7948937 0.2828
6 passive.heating thorax.horn -0.0505 22.86 high 0.7948937 0.0505
6 passive.heating elytra -0.0280 22.86 high 0.7948937 0.0280
6 passive.heating abdomen -0.4571 22.86 high 0.7948937 0.4571 
beetle.passive.mean<-beetle.passive %>%
  group_by(body.part) %>%
  summarise(avg = mean(invert), sd = sd(invert)); beetle.passive.mean
## # A tibble: 4 × 3
##   body.part      avg     sd
##   <chr>        <dbl>  <dbl>
## 1 abdomen     0.224  0.138 
## 2 cephal.horn 0.479  0.405 
## 3 elytra      0.0735 0.0823
## 4 thorax.horn 0.0740 0.0455


6 Control

6.1 Average slope per body part

beet.cont.mean<-beetle.control %>%
  group_by(body.part) %>%
  summarise(avg = mean(invert), sd = sd(invert)); beet.cont.mean
## # A tibble: 4 × 3
##   body.part     avg    sd
##   <chr>       <dbl> <dbl>
## 1 abdomen     0.520 0.439
## 2 cephal.horn 0.745 0.276
## 3 elytra      0.443 0.467
## 4 thorax.horn 0.676 0.306


7 Passive heating

7.1 Does cooling rate differ among body parts?

We are accounting for body mass in the model

lm1<-lm(invert ~ body.part + mass, data = beetle.passive)
check_model(lm1)

checkresiduals(lm1)

## 
##  Breusch-Godfrey test for serial correlation of order up to 8
## 
## data:  Residuals
## LM test = 7.0906, df = 8, p-value = 0.5269
summary(lm1)
## 
## Call:
## lm(formula = invert ~ body.part + mass, data = beetle.passive)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.26053 -0.08271 -0.02570  0.02655  0.65450 
## 
## Coefficients:
##                      Estimate Std. Error t value Pr(>|t|)  
## (Intercept)           0.68076    0.26186   2.600   0.0176 *
## body.partcephal.horn  0.25467    0.11940   2.133   0.0462 *
## body.partelytra      -0.15055    0.11940  -1.261   0.2226  
## body.partthorax.horn -0.15003    0.11940  -1.257   0.2241  
## mass                 -0.02179    0.01183  -1.842   0.0811 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2068 on 19 degrees of freedom
## Multiple R-squared:  0.4968, Adjusted R-squared:  0.3908 
## F-statistic: 4.689 on 4 and 19 DF,  p-value: 0.008396
summary.aov(lm1)
##             Df Sum Sq Mean Sq F value  Pr(>F)   
## body.part    3 0.6570 0.21901   5.121 0.00917 **
## mass         1 0.1452 0.14518   3.395 0.08107 . 
## Residuals   19 0.8126 0.04277                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

There is variation among body parts, with the cephalic horn cooling at a higher rate than other parts. Body mass did not mediate this effect.

tab_model(lm1)
  invert
Predictors Estimates CI p
(Intercept) 0.68 0.13 – 1.23 0.018
body part [cephal.horn] 0.25 0.00 – 0.50 0.046
body part [elytra] -0.15 -0.40 – 0.10 0.223
body part [thorax.horn] -0.15 -0.40 – 0.10 0.224
mass -0.02 -0.05 – 0.00 0.081
Observations 24
R2 / R2 adjusted 0.497 / 0.391 

## 
##  Breusch-Godfrey test for serial correlation of order up to 8
## 
## data:  Residuals
## LM test = 4.9245, df = 8, p-value = 0.7656
## 
## Call:
## lm(formula = invert ~ body.part + mass.cat, data = beetle.passive)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.25640 -0.10392 -0.01390  0.01905  0.73749 
## 
## Coefficients:
##                      Estimate Std. Error t value Pr(>|t|)  
## (Intercept)            0.1869     0.0938   1.993   0.0608 .
## body.partcephal.horn   0.2547     0.1251   2.036   0.0559 .
## body.partelytra       -0.1505     0.1251  -1.204   0.2435  
## body.partthorax.horn  -0.1500     0.1251  -1.200   0.2450  
## mass.catlow            0.1114     0.0938   1.188   0.2496  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2166 on 19 degrees of freedom
## Multiple R-squared:  0.4479, Adjusted R-squared:  0.3316 
## F-statistic: 3.853 on 4 and 19 DF,  p-value: 0.01862
##             Df Sum Sq Mean Sq F value Pr(>F)  
## body.part    3 0.6570 0.21901   4.667 0.0132 *
## mass.cat     1 0.0662 0.06620   1.411 0.2496  
## Residuals   19 0.8916 0.04693                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


8 Plot

Fig1<-ggplot(beetle.passive, aes(x = mass, y = invert, colour=body.part, group=body.part, fill=body.part)) + 
  
  geom_segment(x=12,xend=24,y=beet.cont.mean$avg[1], yend=beet.cont.mean$avg[1], 
               linetype="dashed", color = "#EBC397", lwd=.5, alpha=.45)+
  geom_segment(x=12,xend=24,y=beet.cont.mean$avg[2], yend=beet.cont.mean$avg[2], 
               linetype="dashed", color = "#A6E7FF", lwd=.5, alpha=.45)+
  geom_segment(x=12,xend=24,y=beet.cont.mean$avg[3], yend=beet.cont.mean$avg[3], 
               linetype="dashed", color = "#2C809E", lwd=.5, alpha=.45)+
  geom_segment(x=12,xend=24,y=beet.cont.mean$avg[4], yend=beet.cont.mean$avg[4], 
               linetype="dashed", color = "#9E5F1C", lwd=.5, alpha=.45)+
  
  stat_smooth(formula = "y ~ x", method = "lm", fullrange = T, alpha=0.5,
              show.legend = FALSE, se=F, lwd=.95) +
  
  scale_colour_manual(values=c("#EBC397", "#A6E7FF","#2C809E", "#9E5F1C"))+
  
  geom_point(shape = 21, colour = "black", size = 3, alpha = .9) +
  
  scale_fill_manual(name = "", labels = c("Abdomen","Cephalic horn", "Thoracic horn","Elytra"),
                    values=c("#EBC397", "#A6E7FF","#2C809E", "#9E5F1C"))+
  
  ylab("Cooling rate\n(slope)") +
  xlab("Body mass (g)") +
  
  theme(
    # legend.position=c(0.8,0.85),
    legend.title=element_text(size=11),
    legend.text=element_text(size=10)) +
  
  scale_y_continuous(limits=c(0,1.5), labels=seq(0,1.5,0.5))+
  
  ggtheme(12,0.3)

# Beetle1 <- readPNG(paste0(figuresfolder, "/Beetle inset1.png"), native=T)
# Beetle2 <- readPNG(paste0(figuresfolder, "/Beetle inset2.png"), native=T)
Beetle3 <- readPNG(paste0(figuresfolder, "/Beetle inset3.png"), native=T)

Fig1Inset<- Fig1+
  inset_element(p=Beetle3, 0.55, 0.7,1,.9); Fig1Inset
**Figure 3**. Cooling rate as a function of body mass in *Megasoma gyas*. Body parts are color-coded, highlighting that heat dissipation occurs at different rates among body parts. The solid lines represent the expected relationship between cooling rates and body mass per body part. The dashed lines represent the cooling rates of different body parts measured from control specimens.

Figure 3. Cooling rate as a function of body mass in Megasoma gyas. Body parts are color-coded, highlighting that heat dissipation occurs at different rates among body parts. The solid lines represent the expected relationship between cooling rates and body mass per body part. The dashed lines represent the cooling rates of different body parts measured from control specimens. 

plot(invert ~ as.factor(mass.cat), data = beetle.passive)


9 Validation of results

One individual (ID5) was much lighter than the others in our study. This individual had the cephalic horn with the highest cooling rate in our data set. To understand if the result we obtained was caused by a potential outlier, we will remove ID5 from our data set and re-run our analyses.

beetle.passive.subset <- beetle.passive %>%
  filter(ID != '5' )

9.1 Does cooling rate differ among body parts?

We are accounting for body mass in the model

lm3.1<-lm(invert ~ body.part+ mass, data = beetle.passive.subset)
check_model(lm3.1)

checkresiduals(lm3.1)

## 
##  Breusch-Godfrey test for serial correlation of order up to 8
## 
## data:  Residuals
## LM test = 11.225, df = 8, p-value = 0.1893
summary(lm3.1)
## 
## Call:
## lm(formula = invert ~ body.part + mass, data = beetle.passive.subset)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.16834 -0.04533 -0.01075  0.04340  0.20989 
## 
## Coefficients:
##                       Estimate Std. Error t value Pr(>|t|)   
## (Intercept)           0.123462   0.298289   0.414  0.68480   
## body.partcephal.horn  0.071660   0.062187   1.152  0.26722   
## body.partelytra      -0.174860   0.062187  -2.812  0.01314 * 
## body.partthorax.horn -0.186160   0.062187  -2.994  0.00909 **
## mass                  0.005413   0.013171   0.411  0.68688   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.09833 on 15 degrees of freedom
## Multiple R-squared:  0.6318, Adjusted R-squared:  0.5336 
## F-statistic: 6.434 on 4 and 15 DF,  p-value: 0.003197
summary.aov(lm3.1)
##             Df  Sum Sq Mean Sq F value  Pr(>F)   
## body.part    3 0.24717 0.08239   8.522 0.00153 **
## mass         1 0.00163 0.00163   0.169 0.68688   
## Residuals   15 0.14502 0.00967                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
lm3.2<-lm(invert ~ body.part-1 + mass, data = beetle.passive.subset)
summary.aov(lm3.2)
##           Df Sum Sq Mean Sq F value   Pr(>F)    
## body.part  4 0.8415 0.21037  21.759 4.21e-06 ***
## mass       1 0.0016 0.00163   0.169    0.687    
## Residuals 15 0.1450 0.00967                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Apparently, the cephalic horn and abdomen have cooling rates more similar to the control compared to the elytra and toracic horn. This makes sense, as flight muscles (which are known to play a role in insect thermoregulation) are concentrated in the region of the elytra/thorax.

## 
##  Breusch-Godfrey test for serial correlation of order up to 8
## 
## data:  Residuals
## LM test = 7.8528, df = 8, p-value = 0.448
## 
## Call:
## lm(formula = invert ~ body.part + mass.cat, data = beetle.passive.subset)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.16389 -0.04585 -0.01375  0.04595  0.20931 
## 
## Coefficients:
##                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)           0.24779    0.04546   5.450  6.7e-05 ***
## body.partcephal.horn  0.07166    0.06238   1.149  0.26861    
## body.partelytra      -0.17486    0.06238  -2.803  0.01337 *  
## body.partthorax.horn -0.18616    0.06238  -2.985  0.00926 ** 
## mass.catlow          -0.01535    0.05513  -0.278  0.78449    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.09862 on 15 degrees of freedom
## Multiple R-squared:  0.6295, Adjusted R-squared:  0.5307 
## F-statistic: 6.372 on 4 and 15 DF,  p-value: 0.003335
##             Df  Sum Sq Mean Sq F value  Pr(>F)   
## body.part    3 0.24717 0.08239   8.471 0.00157 **
## mass.cat     1 0.00075 0.00075   0.078 0.78449   
## Residuals   15 0.14590 0.00973                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

9.2 Plot

Fig2<-ggplot(beetle.passive.subset, aes(x = mass, y = invert, colour=body.part, group=body.part, fill=body.part)) + 
  
  geom_segment(x=12,xend=24,y=beet.cont.mean$avg[1], yend=beet.cont.mean$avg[1], 
               linetype="dashed", color = "#EBC397", lwd=.5, alpha=.45)+
  geom_segment(x=12,xend=24,y=beet.cont.mean$avg[2], yend=beet.cont.mean$avg[2], 
               linetype="dashed", color = "#A6E7FF", lwd=.5, alpha=.45)+
  geom_segment(x=12,xend=24,y=beet.cont.mean$avg[3], yend=beet.cont.mean$avg[3], 
               linetype="dashed", color = "#2C809E", lwd=.5, alpha=.45)+
  geom_segment(x=12,xend=24,y=beet.cont.mean$avg[4], yend=beet.cont.mean$avg[4], 
               linetype="dashed", color = "#9E5F1C", lwd=.5, alpha=.45)+
  
  stat_smooth(formula = "y ~ x", method = "lm", fullrange = T, alpha=0.5,
              show.legend = FALSE, se=F, lwd=.95) +
  
  scale_colour_manual(values=c("#EBC397", "#A6E7FF","#2C809E", "#9E5F1C"))+
  
  geom_point(shape = 21, colour = "black", size = 3, alpha = .9) +
  
  scale_fill_manual(name = "", labels = c("Abdomen","Cephalic horn", "Thoracic horn","Elytra"),
                    values=c("#EBC397", "#A6E7FF","#2C809E", "#9E5F1C"))+
  
  ylab("Cooling rate\n(slope)") +
  xlab("Body mass (g)") +
  
  theme(
    # legend.position=c(0.8,0.85),
        legend.title=element_text(size=11),
        legend.text=element_text(size=10)) +
  
  scale_y_continuous(limits=c(0,1.5), labels=seq(0,1.5,0.5))+
  
  ggtheme(12,0.3);Fig2

9.2.1 Grid plot

Notice the different x-axis dimensions.

cowplot::plot_grid(Fig1, Fig2, labels = c('A', 'B'),
          ncol=2, align = "hv")

sessionInfo() %>% pander

R version 4.3.2 (2023-10-31)

Platform: aarch64-apple-darwin20 (64-bit)

locale: en_US.UTF-8||en_US.UTF-8||en_US.UTF-8||C||en_US.UTF-8||en_US.UTF-8

attached base packages: stats, graphics, grDevices, utils, datasets, methods and base

other attached packages: sjPlot(v.2.8.15), patchwork(v.1.2.0), png(v.0.1-8), pander(v.0.6.5), cowplot(v.1.1.2), kableExtra(v.1.3.4), forecast(v.8.21.1), lubridate(v.1.9.3), forcats(v.1.0.0), stringr(v.1.5.1), dplyr(v.1.1.4), purrr(v.1.0.2), readr(v.2.1.5), tidyr(v.1.3.0), tibble(v.3.2.1), ggplot2(v.3.5.1), tidyverse(v.2.0.0), lme4(v.1.1-35.1), Matrix(v.1.6-1.1), performance(v.0.10.8) and scales(v.1.3.0)

loaded via a namespace (and not attached): rlang(v.1.1.3), magrittr(v.2.0.3), tseries(v.0.10-55), compiler(v.4.3.2), mgcv(v.1.9-0), systemfonts(v.1.0.5), vctrs(v.0.6.5), quadprog(v.1.5-8), rvest(v.1.0.3), pkgconfig(v.2.0.3), fastmap(v.1.1.1), backports(v.1.4.1), labeling(v.0.4.3), effectsize(v.0.8.6), utf8(v.1.2.4), rmarkdown(v.2.25), tzdb(v.0.4.0), nloptr(v.2.0.3), xfun(v.0.41), cachem(v.1.0.8), jsonlite(v.1.8.8), highr(v.0.10), sjmisc(v.2.8.9), ggeffects(v.1.3.4), broom(v.1.0.5), parallel(v.4.3.2), R6(v.2.5.1), bslib(v.0.6.1), stringi(v.1.8.3), boot(v.1.3-28.1), lmtest(v.0.9-40), jquerylib(v.0.1.4), estimability(v.1.4.1), Rcpp(v.1.0.12), knitr(v.1.45), modelr(v.0.1.11), zoo(v.1.8-12), parameters(v.0.21.3), splines(v.4.3.2), nnet(v.7.3-19), timechange(v.0.2.0), tidyselect(v.1.2.0), rstudioapi(v.0.15.0), yaml(v.2.3.8), timeDate(v.4032.109), sjlabelled(v.1.2.0), curl(v.5.2.0), lattice(v.0.22-5), quantmod(v.0.4.25), withr(v.2.5.2), bayestestR(v.0.13.1), urca(v.1.3-3), coda(v.0.19-4.1), evaluate(v.0.23), xts(v.0.13.1), xml2(v.1.3.6), pillar(v.1.9.0), insight(v.0.19.7), generics(v.0.1.3), TTR(v.0.24.4), hms(v.1.1.3), munsell(v.0.5.0), minqa(v.1.2.6), xtable(v.1.8-4), glue(v.1.7.0), emmeans(v.1.9.0), tools(v.4.3.2), see(v.0.8.1), webshot(v.0.5.5), mvtnorm(v.1.2-4), grid(v.4.3.2), datawizard(v.0.9.1), colorspace(v.2.1-0), nlme(v.3.1-163), fracdiff(v.1.5-2), cli(v.3.6.2), fansi(v.1.0.6), viridisLite(v.0.4.2), svglite(v.2.1.3), sjstats(v.0.18.2), gtable(v.0.3.4), sass(v.0.4.8), digest(v.0.6.34), ggrepel(v.0.9.5), farver(v.2.1.1), htmltools(v.0.5.7), lifecycle(v.1.0.4), httr(v.1.4.7) and MASS(v.7.3-60)